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Abstract. A simple model is proposed as a possible explanation for the checkerboard pattern of modulations 
in the hole density observed in recent tunneling experiments on underdoped cuprates. Two assumptions 
are made; first, an enhanced hole density near the acceptor dopants and secondly short range correlations 
in the positions of these dopants caused by their electrostatic and anisotropic elastic interactions. Together 
these can lead to a structure factor in qualitative agreement with experiment. 

PACS. 74.25.Jb Electronic structure - 74.72.-h Cuprate superconductors 



Recent studies which measured the tunneling conduc- 
tance into underdoped cuprates found an unexpected spa- 
tial modulation at higher bias voltages with a 4x4 period- 
icity. The most comprehensive studies pQ have been made 
on the Ca2-i E Na ;r Cu02Cl2 compounds where a good sur- 
face can be prepared deep into the underdoped region. 
These scanning tunneling microscopy (STM) studies have 
led to interpretations of the data as signals of superlattices 

' of single holes or hole pairs JZEBGIBJ ■ Since these mate- 
rials are Mott insulators doped randomly with acceptors 
it is not clear if superlattice ordering can persist in the 
presence of the random electric fields of the acceptor ions. 
In this letter we will explore an alternative explanation 
that starts from the assumption of enhanced hole densi- 
ties localized in the vicinity of acceptors and examines the 

' consequences of short range correlations in the positions 
of the acceptors, e.g. the Na-ions in Ca2- x Na K Cu02Cl2. 

Our first assumption of enhanced hole densities close 
to the acceptors is consistent with the standard behavior 
' of doped semiconductors. These undergo a Mott transi- 
tion when the dopant concentration is reduced and hole 
wavefunctions become localized around acceptors. Here we 
make a simple assumption of an enhanced hole density and 
explore the consequences. Then we simulate the positions 
of the Na ions under the combined influence of inter-ion 
Coulomb repulsion and elastic forces, using a Monte Carlo 
algorithm. Together these two effects determine the form 
of the modulations in the hole density and so the STM 
structure factor. As we shall see the calculated structure 
factor shows qualitative agreement with the STM experi- 
ments. 

The outermost layer of the crystal consists of inter- 
penetrating square lattices of CI and (Ca^^Na^) ions. 
The STM tip will couple most effectively to the outermost 



wavefunctions of the surface layer which should be the 3p z 
orbitals of the CI ions. The final state of the injected elec- 
tron (+ve bias) is in the hybridized antibonding band of 



3d x 2 



orbitals centered on Cu sites and the 2p x y O- 



orbitals. Since the CI ions lie directly above the Cu-sites 
there is no matrix element to tunnel into the state with 
d x 2_ y 2 -symmetry centered of the Cu-site directly below. 
There will be however a finite matrix element to couple 
into the 4 nearest neighbor Cu sites generated through the 
overlap of the 3p z -Cl and 2p xy -0 states on neighboring 
Cl-0 pairs. This situation is similar to the case of tun- 
neling through an outermost BiO layer in BSSCO com- 
pounds which has important consequences for the form 
factor of substitutional impurities in the Cu02-plane as 
was pointed out by Martin, Balatsky and Zaanen [7j. Re- 
turning to the case of Ca2-a;Na^Cu02Cl4 the injection of 
an electron (or hole) through the 3p z orbital of a CI ion 
neighboring a Na-ion leads to finite amplitudes for injec- 
tion into the 4 Cu sites which are nearest neighbors (nn) 
to the Cu underneath the CI ion. Of these 4 sites 2 are nn 
sites of the Na ion and 2 are nnn sites (see Fig. 1). 

The next question that arises concerns the tunneling 
probability is whether the tunneling amplitude into these 
4 sites should be added coherently or incoherently. Con- 
sider the case of positive bias or electron injection from the 
tip. Clearly the Coulomb attraction between a hole and 



Fig. 1. Projected view of crystal structure of 
Ca2-iNa^Cu02Cl2 around a doped Na + ion. Black and 
gray balls represent Cu and Ca ions. An electron (or hole) 
injected to CI 2p z -orbital above site is transferred to one of 
the Cu sites 1-4. 
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the Na-acceptor leads to increased occupation of states 
centered on the 4 n.n. Cu sites and to a lesser extent 
on the 8 n.n.n. sites. Explicit form of hole wavefunction 
around a dopant was numerically calculated for the case 
of Li + dopant in ladder system Similar calculation is 
possible for this case but this is beyond the scope of this 
paper. This in turn enhance the hole density on these sites. 
When we calculate the electron spectral function A(k, ui) 
that determines the conductance 



Fig. 3. (a) Fourier transform of Na-ion distribution, (|p(k)| 2 ), 
annealed for the potential V^(R) with ^4=1.0. Ensemble av- 
erage is taken over 10 2 configurations generated by Monte 
Carlo calculations and the results are smoothed by convolut- 
ing with Gaussian function with width 0.037T. The system size 
is 100x100 sites and the Na ion density is 0.10. (b) Fourier 
transform of random initial configurations at the same density 
before annealing. 
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We are interested in an energy range stretching up to sev- 
eral hundred meV. Therefore the relevant matrix element 
will be a incoherent superpositions of tunneling processes 
into the 4 Cu sites which are n.n. to the CI site. Here <JQ h ) 
is the one- hole groundstate with energy E}^ and l^" 11 } 
are eigenstates of a stoichiometric system, and /i(x, x') is 
the amplitude of the hole wave function centered on the 
Cu-site x', at the CI site at position, x. For calculating 
the spectral function, we average the energy over the al- 
lowed range for the injected electron so that we may treat 
the sum over the one-hole eigenstates as a sum over a 
complete set of states. Under these incoherent conditions 
the (A(x, w)) w _ av ~ J2 X > /i 2 (x, x')n(x') where n(x') is the 
hole density at x'. In the neighborhood of a Na + -ion the 
Coulomb attraction leads to an enhanced hole density on 
the inner ring of 4 nearest neighbor sites and to a lesser 
extent on the outer ring of 8 next nearest neighbor sites 
(see Fig. 1). This in turn will lead to an enhanced tunnel- 
ing signal on the inner ring (/i) and the outer ring (j^)- 
The Fourier transformed cluster form factor is then 



/(k) = 4/i cos y cos y 



+4/ 2 



cos cos — 

2 2 



COS COS - 

2 2 



(2) 



We do not attempt detailed calculations of /i , and fi but 
treat them as phenomenological parameters. Note in this 
letter we will concentrate on electron injection only and 
so our calculations are for A(lu), where to is an energy of 
order 100 [meV] or above, where the experimental data 
do not depend on the energy so much. Strictly speaking 
this should be determined from experiment by normalizing 
the spectra in the -ve bias (hole injection) side rather than 
the +ve bias side chosen by Hanaguri et al.[7|. An explicit 
proposal along these lines has been made by Randeria and 
coworkers [H]- 



Fig. 2. Typical configuration of Na ion distribution simulated 
by Monte Carlo calculations. The system size is 100x100 sites 
and the Na ion density is 0.10. Strength of elastic component 
of Na-Na interaction is (a) A=0.5, and (b) A=1.0. 



The STM signal intensity is given by 

J(x) = J d 2 R P (R)/(x-R) (3) 

where /(x — R) is the intensity of a cluster at CI site x 
around a Na-acceptor at R and it is the inverse Fourier 
transform of /(k) in Eq. J5J. p(R) is the distribution 
function for the Na-ions. At the relevant densities x ks 
0.1 the concentration of the Na ions in the outermost 
Ca2- x Na x Cu02Cl2 layer is not in the dilute regime so 
that it is possible that the Na-Na interactions lead to 
short range order in the positions of the Na acceptors. 
To estimate the effects of these interactions on the distri- 
bution we simulated the annealing of an initial random 
2-dimensional distributions through inter-ion interactions 
using a Monte Carlo algorithm. The interaction has a 
Coulomb repulsion component and an elastic component 



(4) 



The former is radially symmetric but the second term is in 
general anisotropic. We estimate it following EshelbyfTU] 
for a 2-dimensional square lattice and we ignore effects 
associated with the fact that the Na-acceptors lie both 
above and below the topmost CuO-planc. Here 9 is the 
angle of vector R measured from the principal axis of the 
lattice. The intensity A — 2ttc 2 i](X + /i) /(A + 2/i) is given 
by Lame constants, A = X X x,yy and fi = \ X y,xy, and the 
tetragonal anisotropy, rj = X xx ,xx — (A + 2/i), with c being 
a constant proportional to local volume change due to im- 
purity doping. This is normalized such that the coefficient 
of the Coulomb term is set to be unity. 

We have investigated stable configuration of doped 
Na + -ions with the above interaction potential. For this 
simulated annealing, we have used a Monte Carlo algo- 
rithm. Starting from a random initial configuration, we 
choose a Na-ion in turn and find a vacant site randomly. 
We then calculate how much the interaction energy changes, 
AE, if the chosen ion is to move. If AE < 0, we al- 
ways move that ion to achieve a stabler configuration, 
whereas even if AE > 0, we accept the new configura- 
tion with the Boltzmann probability er AE l T at temper- 
ature T. This update process is applied to all the ions in 
turn, and the procedure is repeated many times to obtain 
a stable ion configuration. Typically, we have set the tem- 
perature T — 0.01 — 0.05, and repeated the Monte Carlo 
procedure 10 3 — 10 4 steps for each initial configuration. We 
have checked that obtained results do not depend on these 
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Fig. 4. Structure factor of holes localizing around doped Na 
ions, (|J(k)| 2 ). Ensemble average is taken over 10 2 confi gura- 
tions generated by Monte Carlo calculations and the results are 
smoothed by convoluting with Gaussian function with width 
0.037T. The system size is 100 x 100 sites and the Na ion density 
is 0.10. Coupling constant is (a) A=0.5, (b) A=1.0, and (c) 
A=0.0. 

parameters sensitively. In Fig. 2 we show typical distri- 
butions resulting from Monte Carlo annealing of samples 
with 10 4 sites and 10 3 ions. Periodic boundary conditions 
are imposed to both directions and the interaction V(R) 
is tailored for the boundary conditions by using the mirror 
charge technique. 

Once the annealed configuration is obtained, we calcu- 
late the structure factor of Na ions by taking the Fourier 
transform. Result for 10% doping and A=1.0 is plotted 
in Fig. 3(a), and the Fourier transform of random initial 
configurations is also shown in Fig. 3(b) in comparison. 
We have also checked that ensemble average over a larger 
set of samples (e.g., 10 3 ) gives essentially the same results. 
The Na-ion structure factor is averaged over initial con- 
figurations, and then smoothed by Gaussian broadening 
with width 0.037T. The final structure factor is calculated 
by multiplying it with the form factor of hole density lo- 
calizing around Na ions 

|/(k)| 2 = |p(k)| 2 |/(k)| 2 . (5) 



vicinity of an acceptor which hopefully can be explicitly 
tested in experiment. 
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In Fig. 4 (a) and (b), we show representative sets of data 
for two values of the relative strength of the elastic and 
Coulomb interactions. The data for the case of pure Coulomb 
interaction are also shown in Fig. 4(c) in comparison. 
We see that the main characteristic of the checkerboard 
pattern observed in the experiments, [J] namely peaks at 
wavevectors (n/2,0) and (0,7r/2) are nicely reproduced 
when the anisotropic elastic energy is considered A > 0, 
especially for the case with stronger elastic interactions. 
The stronger tendency for the Na-ions to line up along the 
crystal axes in these simulations is reflected in the struc- 
ture factor by the suppression of weight along the diagonal 
leading to the peaks being concentrated near the vertical 
and horizontal axes. 

Our calculations show that at least the broad features 
of the data can be reproduced by our two assumptions, en- 
hanced hole density around the acceptor and short range 
correlations in the positions of the Na- acceptors. Such 
short range order necessarily induces a corresponding mod- 
ulation in the hole distribution. In addition the hole den- 
sity will be distributed over a cluster of sites around each 
acceptor. Both effects influence the final structure factor. 
The STM data contains a lot of more detailed spectro- 
scopic information, e.g. small peaks in the conductance. 
These features are beyond the scope of the calculations 
presented here. The STM experiments give us really de- 
tailed information on an unprecedented scale but for now 
at least, we can only analyze the main features. A key 
point in our model is the enhanced hole density in the 
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